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Positioning  of  Large  Surface  Vessels  using  Multiple  Tugboats 

D.  Braganza^,  M.  Feemster,  and  D.  Dawson 


Abstract —  In  this  paper,  the  positioning  of  large  surface  ves¬ 
sels  using  multiple,  autonomous  tugboats  is  investigated.  Specif¬ 
ically,  the  paper  investigates  the  development  of  two  control 
strategies:  i)  an  adaptive  position  controller  that  compensates 
for  select  system  parameters  such  as  mass  and  drag  coefficients 
and  also  for  an  unknown  thrust  configuration  matrix  that 
is  the  result  of  unknown  tugboat  locations  and  ii)  a  robust 
position  control  strategy  that  does  not  require  transmission 
of  information  between  the  tugs.  The  control  design  for  ii)  is 
facilitated  by  strategic  placement  of  the  tugs  about  the  vessel 
hull  such  that  the  resulting  thrust  matrix  is  upper  triangular. 
Performance  of  both  controllers  is  investigated  by  simulations 
for  a  scaled  model  of  a  surface  vessel. 

I.  Introduction 

THE  positioning  of  large  surface  vessels  such  as  barges, 
platforms,  and  ships  during  low  speed  maneuvers  in 
such  confined  environments  as  harbors  is  a  daunting  task. 
Due  to  the  low  speed  operation,  large  vessels  that  are  steered 
by  rudders  typically  cannot  generate  sufficient  control  influ¬ 
ence  over  their  position  and  orientation  for  delicate  docking 
maneuvers;  therefore,  multiple  tugboats  are  required  to  assist 
in  positioning  these  vessels.  In  this  scenario,  tugboats  are 
strategically  positioned  along  the  hull  of  the  surface  vessel 
and  their  efforts  are  coordinated  such  that  their  thrust  inputs 
maneuver  the  vessel  along  a  desired  path.  This  manual 
manipulation  of  the  surface  vessel  with  human  operators 
within  the  process  requires  precise  coordination  via  radio 
communication  between  all  involved  tugboats.  Therefore,  its 
performance  is  dependent  on  the  skill  of  the  tugboat  opera¬ 
tors  as  well  as  the  reliability  of  the  communication  network. 
With  the  advent  of  new  propulsion  technologies  such  as 
tunnel  thrusters  and  accurate  differential  global  positioning 
system  (DGPS),  operators  now  posses  new  tools  to  aid  in  this 
task;  however,  complete  automation  coupled  with  reduced 
reliance  on  the  communication  network  between  the  tugs  of 
this  process  poses  a  challenging  problem. 

In  recent  years,  researchers  have  been  focusing  on  various 
marine  control  problems  for  surface  vessels  [  1  ]  which  can  be 
loosely  classified  in  the  following  areas:  i)  position  tracking 
(this  includes  trajectory  tracking  and  path  following  also 
known  as  maneuvering  [2],  [3]),  ii)  way  point  tracking 
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[4],  iii)  development  of  autopilots  for  steering  [5],  iv)  dy¬ 
namic  positioning  or  station  keeping  [6],  [7],  and  recently 
v)  formation  control  [8],  [9].  The  problem  of  uncertainty 
in  the  dynamic  model  has  also  been  addressed  with  the 
development  of  adaptive  controllers  (for  example  [3],  [10]) 
for  systems  where  the  structure  of  the  dynamic  model  is 
known  and  fuzzy/neural  network  control  architectures  [5], 
[11],  [12]  for  systems  whose  dynamic  model  structure  is 
not  known.  Although  the  positioning  problem  for  fully 
actuated  surface  vessels  has  been  considered  in  [6],  [7],  these 
works  typically  assume  that  the  thrust  configuration  matrix  is 
completely  known.  In  [13],  station  keeping  experiments  were 
conducted  for  an  over-actuated  ship  model,  and  a  constrained 
optimization  algorithm  was  used  to  allocate  the  control  to  the 
thrusters  subject  to  saturation  limits  and  other  constraints. 
However,  this  approach  also  requires  the  thrust  configuration 
matrix  be  known  a  priori. 

Based  on  the  literature  review,  it  appears  that  two  scenarios 
within  the  manipulation  of  surface  vessels  with  autonomous 
tugboats  may  be  further  matured:  i)  compensation  for  an 
unknown  thrust  configuration  matrix  (i.e.,  the  locations  of 
the  tugs  around  the  vessel  hull  are  not  known),  and  ii)  a 
control  strategy  which  requires  no  communication  between 
the  tugs.  In  regards  to  i),  many  current  surface  vessel  control 
strategies  do  not  impose  the  additional  design  constraint  of 
an  unknown  thrust  configuration  matrix.  This  omission  may 
be  attributed  to  the  fact  that  previous  ship  controllers  were 
designed  assuming  that  the  vessel  under  consideration  has  the 
ability  to  generate  its  own  thrust  from  known  actuators  posi¬ 
tioned  at  known  locations  and  not  from  externally  connected 
tugboats.  As  a  result  of  this  unknown  thrust  matrix,  many 
commonly  utilized  force  allocation  (commutation)  strategies 
may  be  rendered  ineffective  to  distribute  the  required  thrust 
command  among  the  connected  tugs.  Therefore,  an  adaptive, 
asymptotic  position  tracking  controller  that  compensates  for 
an  unknown  thrust  configuration  matrix  as  well  as  such  select 
system  parameters  as  vessel  mass  and  drag  coefficients  is 
developed. 

In  regards  to  ii),  reduced  communication  between  the  tug¬ 
boats  holds  particular  interest  in  naval  warfare  applications 
such  as  utilizing  a  swarm  of  autonomous  tugs  to  recover  a 
ship  from  a  hostile  environment.  Since  each  tugboat  operates 
in  an  independent  manner,  communication  between  the  tugs 
is  not  susceptible  to  potential  electronic  jamming  which 
may  render  other  coordinated  force  allocation  strategies 
ineffective  .  In  order  to  address  this  omission  of  information 
between  the  tugs,  the  design  for  this  case  is  facilitated  by 
strategic  placement  of  the  tugs  around  the  vessel  hull  such 
that  the  thrust  configuration  matrix  is  in  an  upper  triangular 


form  thereby  allowing  for  a  hierarchical  design  technique. 
In  essence,  “upper  level”  tugs  the  treat  the  thrusts  from  the 
"lower  level”  tugs  as  disturbances  which  are  duly  compen¬ 
sated  for  via  a  neural  network  feedforward  component  [14] 
and  a  robust  control  structure. 

In  addition  to  the  above  constraints,  the  control  strategies 
for  both  cases  are  designed  in  a  manner  to  accommodate 
manipulation  of  large  vessels  where  mooring  of  the  tugboats 
is  not  tractable  (i.e.  secure  attachment  of  the  tugs  to  the  hull 
cannot  be  achieved).  That  is,  the  control  algorithm  for  each 
tug  is  designed  such  that  only  a  positive  or  “pushing”  thrust 
is  requested  (i.e.,  the  tugs  cannot  pull  on  the  vessel).  As  a 
result,  a  minimum  of  six  tugboats  are  required  for  the  three 
degree  of  freedom  positioning  problem  (x,  y  positioning  and 
if)  heading  control).  At  this  stage  it  should  be  mentioned  that 
one  major  underlying  assumption  of  the  control  development 
is  the  static  positioning  of  the  tugboats.  That  is,  it  is  assumed 
that  once  a  tugboat  has  made  contact  with  the  hull  it  remains 
in  contact  and  does  not  slip  along  or  change  incident  angles 
on  the  hull.  Of  course,  this  assumption  may  be  viewed  by 
some  as  not  practical  or  too  restrictive;  however,  this  assump¬ 
tion  allows  focus  on  the  adaptive  and  reduced  communication 
aspects  of  the  control  problem.  As  seen  in  the  subsequent 
sections,  the  tugboat  dynamics  have  not  been  included  within 
the  scope  of  this  paper.  As  a  result,  future  research  work  will 
target  the  inclusion  of  the  dynamics  of  the  tugboats  in  an 
effort  to  remove  this  assumption. 

The  paper  is  structured  in  the  following  manner.  In  Sec¬ 
tion  II,  the  thrust  configuration  matrix  is  studied  through 
a  force  decomposition  analysis  and  tugboat  placement  is 
selected.  The  3  DOF  kinematic  and  dynamic  model  of  a 
vessel  influenced  by  a  swarm  of  tugboats  is  presented  in 
Section  III.  In  Section  IV,  an  adaptive  controller  is  designed 
to  compensate  for  an  unknown  thrust  configuration  matrix 
under  the  assumption  of  full  communication  between  the 
tugboats.  Subsequently  in  Section  V,  a  robust  positioning 
controller  is  developed  for  the  case  when  communication 
between  the  tugboats  is  not  available.  Simulation  results 
for  both  controllers  are  presented  in  Section  VI  to  verify 
their  performance.  Conclusions  and  future  research  efforts 
are  discussed  in  Section  VII. 


where  Ui(t)  =  [ula,ulb,ii2a,U2b,u3a,U3b}T  £  R6  denotes 
a  unidirectional  thrust  input  vector  from  the  six  tugboats  (i.e., 
ui  s>  0)  and  B\  £  R3x6  represents  the  thrust  configuration 
matrix  defined  as  follows 
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where  Cg  =  cos(6),  and  Sg  =  sin(0).  In  order  to  account 
for  the  tugboat’s  unidirectional  influence  on  the  barge,  the 
tugboats  are  arranged  (i.e.  placed)  in  an  opposing  pair  fashion 
about  the  vessel  hull  as  follows 

0^2 b  —  0-2 a  T"  7T>  Lx2b  —  L/x2 a,  Ly2b  Ly2a 

b  —  £^3  a  7T ,  Lx  3b  —  Lx  3a,  Ly3b  —  Ly3a 

which  then  results  in  the  following  force  expression 


II.  Force  Decomposition  and  Commutation 
Strategy 

The  vessel  of  Figure  1  is  assumed  to  be  an  unactuated  (i.e. 
a  barge  or  disabled  ship)  whereby  the  only  thrust  affecting 
the  vessel  is  generated  from  the  six  tugboats  in  contact  with 
the  vessel’s  hull.  Each  tugboat  is  only  allowed  to  apply  a 
unidirectional  thrust  input  to  the  vessel  (i.e.,  “push”  on  the 
hull).  In  addition,  it  is  assumed  that  the  position  and  incident 
angle  of  the  tugboats  are  constant  (the  tugboat  dynamics  and 
contact  issues  such  as  slipping  are  not  considered  within  the 
scope  of  this  paper).  The  total  forces  and  moments  acting  on 
the  vessels  center  of  mass  (COM)  F  =  [Fx,  Fy,  MZ]J  £  R3 
is  the  result  of  the  combined  efforts  of  the  six  tugboats  and 
are  related  by  the  following  expression 

(1) 


F  =  BU 


(4) 


where  the  U(t)  =  [ui,U2,u3]r  £  R3  represents  the  bidi¬ 
rectional  thrust  vector,  Ui  =  Uia  —  Uib,Wi  =  1,2,3  and  the 
resulting  thrust  configuration  matrix  B  £  R3x3  is  specified 
as 


B  = 


1 


Ca2a 
Ca  3a 


0 


$a2a 

$a3a 


0 

IJy2aCa2a  ~  IJx2aSa2a 

—  Ly3aCa3a  IJx3aSa3a 


T 


(5) 


Remark  1:  Upon  closer  inspection  of  the  structure  of  (5), 
the  angle  of  incidence  of  tugboat  U2a  can  be  selected  as 
follows  such  that  the  thrust  configuration  matrix  B  is  placed 
in  upper  triangular  form 


0-2a 


tUTl  (Ly2a/ Lx2a)  • 


F  =  B\  U i 


(6) 


Being  able  to  place  B  into  this  form  allows  for  a  robust  con¬ 
trol  design  to  be  developed  that  requires  no  communication 
between  each  tug.  In  essence,  “upper  level”  tug  pairs  such  as 
u\  view  the  thrust  inputs  from  the  “lower  level”  tugs  of  a  2 
and  U3  as  disturbances  within  the  x  translational  dynamics. 

Remark  2:  It  is  interesting  to  note  that,  the  locations  of  the 
tugboats  specified  in  (3),  can  be  specified  a  priori  such  that, 
no  communication  is  required  between  the  tugboats  while 
maneuvering  the  vessel.  Also,  each  tugboat  will  only  require 
its  own  position  and  orientation  relative  to  the  vessels  COM 
be  known. 

In  the  subsequent  sections,  a  bidirectional  thrust  input 
U(t)  will  be  designed  which  promotes  the  position  tracking 
objective  under  each  case’s  design  constraints.  Once  a  U(t ) 
is  defined,  the  following  commutation  strategy  [15],  can  be 
utilized  to  specify  each  tug’s  unidirectional  thrust  command 


where  i  =  1,2,3  and  £  R+  represents  a  control  gain 
selected  to  ensure  that  the  thrust  inputs  Ujj  (t)  are  non-zero, 
thus  preventing  the  tugboats  from  losing  contact  with  the 
vessel. 

III.  System  Models 

The  kinematic  and  dynamic  equations  for  a  3  DOF  un¬ 
actuated  vessel  manipulated  by  six  external  tugboats  can  be 
written  as  follows  [  1  ] 

77  =  R(ip)u  (8) 

M(v)v  +  D(v)  =  BU  (9) 

where  M  (v)  £  R3x3  represents  the  mass/inertia  matrix  that 
includes  such  effects  as  added  mass,  D(v)  £  R3  represents 
the  drag,  damping,  and  other  parasitic  effects,  B  £  R3x3 
and  U(t)  £  R3  denote  the  thrust  configuration  matrix  and 
thrust  input  vector  of  (4),  77(f)  =  [x(t),y(t),ip(t))T  £  R3 
represents  the  composite  inertial  position  (x,y)  and  heading 
(t/0  vector  of  the  vessel,  is(t)  =  [ u ,  v,  r]T  £  R3  represents  the 
body  fixed  velocity  signals,  and  R(ip)  €=  SO (3)  denotes  the 
rotation  matrix  that  translates  vessel  coordinates  into  inertial 
coordinates  and  is  defined  as  follows 

cos(ip)  —sin(i/})  0 

R(ip)  =  sin\ip)  cosipfi)  0  .  (10) 

0  0  1 

The  subsequent  development  exploits  the  following  proper¬ 
ties  of  the  system  model  presented  in  (8)-(9) 

Property  1:  The  mass/inertia  matrix  M(y)  and  the  drag 
vector  D(v)  are  bounded  if  v{t)  is  bounded. 

Property  2:  The  mass/inertia  matrix  M(y)  is  positive- 
definite  and  symmetric. 

Property  3:  The  rotation  matrix  lihp)  satisfies  the  follow¬ 
ing  properties 

RTy>m)  =  h,  IW)II  =  1,W>  nn 

R(ip)  =  ~'ipSRT(i/t)  1  ; 


where  the  matrix  skew-symmetric  matrix  S  £  R3x3  is  given 
by 

0  1  0 

-10  0.  (12) 
0  0  0 

In  addition  to  the  above  system  properties,  the  subsequent 
control  strategies  are  designed  under  the  following  assump¬ 
tions 

Assumption  1:  Each  tugboat  knows  its  own  position  Li  = 
(. Lxi,Lyi )  and  orientation  a* Vi  =  la,  lb,  ..,36,  with  respect 
to  the  vessel  center  of  mass. 

Assumption  2:  Each  tugboat  has  access  to  measurements 
of  inertial  position  77(f)  and  inertial  velocity  77(f)  of  the 
vessel  center  of  mass.  This  assumption  can  be  supported  in 
part  by  Assumption  1.  That  is,  each  tugboat  is  assumed  to 
have  knowledge  of  its  position  and  orientation  with  respect 
to  the  vessels  COM.  Since  each  tugboat  also  has  access 
to  measurements  of  its  inertial  position  (typically  available 
through  DGPS),  then  using  the  geometric  relationship  of  its 
hull  location  relative  to  the  vessels  center  of  mass,  each 
tugboat  can  formulate  the  inertial  position  and  orientation 
of  the  vessel. 

Assumption  3:  The  tugboats  remain  in  contact  with  the 
vessel’s  hull  at  all  times.  They  also  remain  at  fixed  locations 
with  respect  to  vessel’s  center  of  mass,  and  do  not  change 
incident  angles  with  respect  to  the  hull  thereby  resulting  in 
a  constant  thrust  configuration  matrix  B. 

IV.  Case  I:  Adaptive  Control  Design 

In  this  section,  an  adaptive  position  control  strategy  is 
developed  for  the  vessel  dynamics  of  (9)  with  the  main 
objective  of  promoting  vessel  tracking  along  a  desired  path 
77 d(t)  £  R3  (it  is  assumed  that  the  selected  vessel  trajectory 
is  sufficiently  smooth  in  the  sense  that  77 d(t),  i/,/(f), 

77 d  (f)  £  Coo).  In  addition,  the  design  must  compensate  for 
the  fact  that  the  tugboats  are  dispersed  at  unknown  locations 
around  the  hull.  As  a  result,  the  thrust  configuration  matrix 
of  (5)  is  partially  unknown  (one  must  remember  that  each  tug 
is  assumed  to  know  its  location  with  respect  to  the  vessel’s 
center  of  mass).  In  addition,  the  control  design  will  also  com¬ 
pensate  for  the  coefficients  associated  with  the  mass  matrix 
M(y)  and  the  parasitic  effects  encapsulated  in  D(v).  As  a 
result  of  the  uncertainty  in  B ,  the  coordinated  control  input 
U{t)  consists  of  coupling  terms  of  all  the  tugs;  therefore, 
communication  between  the  tugboats  is  required  ( i.  e. ,  this 
implies  the  calculation  of  the  individual  tugboat  thrust  inputs 
Ui{t)  will  be  in  some  part  dependent  on  information  from 
other  tugboats). 

A.  Transformed  System  Model 

In  order  to  ease  the  development  of  the  adaptive  controller, 
the  system  model  presented  in  (8)-(9)  is  rewritten  in  the 
following  form  [1] 

M*  {rj)ri  +  D*(v,  rj,  rj)  =  JU  (13) 

where  the  transformations  v  =  f?T(  tp)f]  and  1 >  =  RT — 
if) S Rt  (7/0?)  were  utilized,  The  transformed  system  matrices 


M*(rj)  G  R3x3,  D*(y,  r],  ?))  G  R3,  J(rj)  G  R3x3  in  (13)  are 
defined  in  the  following  manner 

M*  =  R(^)M(u)Rt  (V>) 

D*  =  R(ip)D(v)  -  i’R(i/j)M(v)SRT(i/j)i)  (14) 

It  should  be  noted  here  that  the  positive  definiteness  and 
symmetric  attributes  of  the  original  mass/inertia  matrix  M  (y) 
are  preserved  with  M*(j])  in  (14)  (i.e.,  M*(rj)  is  positive 
definite  and  symmetric).  In  addition,  it  is  assumed  that 
parametric  uncertainty  associated  with  the  system  matrices 
M*(ji)  and  D*(v,  77,  ?))  is  linearly  parameterizable.  Also,  it 
is  assumed  that  the  control  input  term  J(rj,di)U  can  be 
manipulated  into  the  following  form 

JU  =  Y\6\  VI7  Gt3  (15) 

where  kj  (U.  q)  G  R3xpi  represents  a  known,  measurable 
regression  vector  and  9 1  G  RPl  represents  a  vector  of 
unknown,  positive  constants.  In  order  to  prevent  potential 
control  singularities,  it  is  assumed  that  the  upper  and  lower 
bounds  for  the  parameter  9 1  are  known  such  that  J{r),9 1) 
is  invertible.  The  bounds  for  each  element  in  9\  can  be 
calculated  as  follows 

£u  <  Ou  <  On  (16) 

where  Q.u,0u  G  R  denote  the  known  lower  and  upper 
bounds,  respectively,  for  the  ith  element  of  9\ . 


where  a  G  R+  represents  a  constant  control  gain.  After 
taking  the  time  derivative  of  r(t),  multiplying  through  by  the 
mass  matrix  and  substituting  in  the  system  dynamics 

of  (13),  the  open-loop  dynamics  for  r(t)  can  be  expressed 
by  the  following  expression 

M*f  =  -^M*r  +  Y202  +  Yi91  +  JU  (19) 

where  the  definition  of  (15)  has  been  utilized,  9\  =  9\  — 
9\  G  Rp3  represents  the  parameter  estimation  error  with 
9 1  it)  denoting  the  vector  of  parameter  estimates,  the  terms 
7}M*(ri)r(t)  and  </(??,  9i)U{t)  have  been  added/sub  traded 
to  the  right  hand  side  of  (19),  Y262  G  R3  captures  the 
uncertainties  associated  with  the  elements  of  M*  ( rj )  and  D* 
and  is  explicitly  defined  as 

3-202  =  \M*r  +  M*  (ae  -  r/d)  -  D*  (20) 

where  G  R3xp2  is  a  known,  measurable  regression 

vector  and  $2  G  RP2  is  a  vector  of  P2  unknown  constants. 
J(r],9 1)  G  R3x3  represents  a  dynamic  estimate  of  the  input 
related  matrix  J{rj,9 1)  and  is  constructed  from  the  following 
on-line  update  law 

0i  =  proj  {y}  (21) 

where  the  auxiliary  term  y(t)  G  RPl  is  defined  as 

y  =  T1Y1Tr  (22) 


B.  Tracking  Error  System  Development 

To  quantify  the  positioning  objective,  a  position  tracking 
error  signal  e{t)  G  R3  is  defined  as  follows 

e  =  V~Vd  (17) 

where  the  desired  trajectory  r]d(t)  and  its  first  two  derivatives 
rld(^)i  Vd(t)>  Vd(t)  ^  are  known  by  all  the  tugboats  and 
are  bounded  functions  of  time. 

Remark  3:  Since  rjd(t)  is  explicitly  dependent  on  time,  the 
tugboat  group  must  be  synchronized  in  some  fashion.  The 
most  simplistic  method  is  to  reset  clock  counters  on  all  the 
tugs  upon  deployment;  however,  one  would  need  to  know 
when  all  tugs  have  established  contact  in  order  to  initiate 
the  maneuvering  control  strategy.  Current  research  efforts 
are  underway  that  focus  on  the  approach  problem  of  the 
autonomous  tugboats.  Since  some  tugs  may  be  required  to 
circumnavigate  around  to  the  far  side  of  the  ship,  contact  time 
may  vary  among  group  members;  therefore,  one  portion  of 
the  focus  of  this  research  is  the  estimation  of  contact  time 
of  all  tugs  from  known  initial  coordinates.  One  may  con¬ 
sider  broadcasting  a  global  initiate  command;  however,  it  is 
desired  to  keep  wireless/radio  communication  at  a  minimum 
so  as  to  reduce  susceptibility  to  electronic  interference. 

In  order  to  ease  the  design  and  analysis  of  the  second 
order  system,  the  filtered  tracking  error  signal  r(t)  G  R3  is 
defined  as  follows 

r  =  e  +  ae  (18) 


with  Ti  G  RPl  x  1  denoting  a  constant  positive  diagonal  gain 
matrix.  The  function  proj{y}  introduced  in  (21)  is  defined 
in  the  following  manner 


proj{yi}  =  < 
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(23) 


(24) 


where  yi(t)  G  M  denotes  the  ith  component  of  y(t),  and 
9u(t)  G  M  denotes  the  ith  component  of  9\{t).  The  above 
projection  algorithm  is  utilized  to  guarantee  that  </(??,  0i) 
is  invertible  by  ensuring  the  parameter  estimates  remain 
bounded  in  the  sense  that  9_u  <  0u{t)  <  9 u.  For  further 
details  of  the  projection  based  update  algorithm  the  reader 
is  referred  [16]). 

Remark  4:  At  this  point,  it  should  be  noted  that  the 
projection  based  update  law  of  (23)  is  the  source  for  the  re¬ 
quirement  of  communication  between  tugs.  Since  placement 
of  tugboats  to  ensure  that  </(??,  9\)  is  diagonal  is  difficult  to 
achieve,  the  regression  vector  Y\(U,rf)  of  (15)  is  comprised 
of  combinations  of  the  tugboat  thrusts  thereby  requiring  the 
passing  of  their  thrust  magnitude  among  all  the  other  tug 
members. 


Based  on  (19)  and  the  subsequent  stability  analysis  the 
thrust  input  U(t)  is  designed  in  the  following  manner 

U  =  J-1  (-Y262  -  Kr )  (25) 

where  K  £  R3x3  denotes  a  positive  constant  diagonal  gain 
matrix  and  d2(t)  £  RP2  denotes  the  parameter  estimate  for 
02  and  is  calculated  from  the  following  dynamics 

d2  =  T2F2t  r.  (26) 

where  £  RP2XP2  represents  a  constant  positive  diagonal 
gain  matrix.  After  substituting  (25)  into  (19),  the  closed  loop 
tracking  error  dynamics  for  r(t)  can  be  expressed  as 

M*f  =  -Kr  -  i M*r  +  Y^i  +  Y292  (27) 

where  92{t)  =  02  —  02  £  RP2  denotes  the  parameter 
estimation  error. 

C.  Stability  Analysis 

Theorem  1:  Given  the  dynamic  model  of  (13),  the  thrust 
control  input  U(t)  specified  in  (25)  coupled  with  the  pa¬ 
rameter  update  laws  of  (21)  and  (26)  ensure  that  the  vessel 
tracking  error  signal  e(t)  is  driven  asymptotically  to  zero  in 
the  sense  that  ||e(i)||  — >  0  as  t  — >  oo. 

Proof:  See  Appendix  I.  ■ 

V.  Case  II:  No  Communication  Between  Tugboats 

In  order  to  eliminate  the  passing  of  information  among 
tugboats  which  is  required  in  the  adaptive  control  design,  a 
robust  control  strategy  is  developed  for  the  position  tracking 
objective  that  does  not  require  the  presence  of  a  communi¬ 
cation  network.  However  in  order  to  achieve  this  goal,  the 
tugboats  are  required  to  be  strategically  placed  around  the 
hull  in  the  manner  of  (6)  such  that  the  thrust  configuration 
matrix  B  of  (5)  is  upper  triangular.  If  this  tug  placement 
is  achieved  and  maintained,  then  a  hierarchical  approach  is 
followed  whereby  “upper  level”  tugs  treat  thrust  inputs  from 
“lower  level”  tugs  as  disturbances  (this  of  course  motivates 
the  utilization  of  a  robust  control  structure).  As  in  the 
adaptive  control  development,  each  tug  is  aware  of  its  own 
location  with  respect  to  the  vessel’s  center  of  mass,  is  pre¬ 
programmed  with  the  vessel’s  desired  maneuver  trajectory 
r] d(t),  and  has  the  ability  to  calculate  the  vessel’s  inertial 
position/velocity  signal  77(f)  and  77(f),  respectively. 

A.  Error  System  Development 

To  promote  the  control  objective,  the  following  position 
tracking  error  signal  z\(t)  £  R3  is  defined  as 

21  =  RT  (V  -  Vd)  (28) 

where  R(tjj)  £  R3x3  is  the  rotation  matrix  previously 
defined  in  (10),  and  the  desired  trajectory  and  its  first  three 
derivatives  r]d(t),  Vd(t),  p  d(t)  £  R3  are  known  by  all 
the  tugboats  and  are  bounded  functions  of  time.  In  order 
to  ease  the  design  and  analysis,  the  filtered  error  signal 
z2(t)  £  R3  is  defined  in  the  following  manner 

22  =  Rt  (77  -  f]d)  +  Kfzi  (29) 


where  K\  £  R3x3  denotes  a  positive  constant  diagonal 
gain  matrix.  From  (28),  the  open-loop  position  tracking  error 
dynamics  for  z\(t)  can  be  formulated  as 

z\  —  — K\Zi  —  tpSzi  H-  z2  (30) 

where  S'  is  a  skew-symmetric  matrix  previously  defined  in 
(12),  and  the  expressions  of  (1 1)  and  (28)  have  been  utilized. 
After  taking  the  time  derivative  of  (29)  and  multiplying  both 
sides  of  the  equation  by  M{y),  the  following  expression  can 
be  obtained  for  the  open-loop  filtered  tracking  tracking  error 
dynamics 

Mz2  =  Mb  +  Mjt{  Ki  Zl  -  RTild  }  (31) 

where  the  kinematic  relationship  of  (8)  has  been  employed. 
After  adding  and  subtracting  the  terms  Zi(t),  \ M{y)z2{t )  to 
(31),  and  substituting  in  (9),  the  following  open  loop  error 
system  for  z2(t)  is  obtained 

Mz2  =  -  i Mz2  -  zx  +  N  +  BU  (32) 

where  the  auxiliary  term  N(-)  £  R3  is  defined  in  the 
following  manner 

N  =  X-Mz2  +  Zl-D  +  Mjf  {KlZl  -  RTild }  .  (33) 

Based  on  (32)  and  the  subsequent  stability  analysis,  the  thrust 
input  vector  U (t)  is  designed  as 

U  =  —B  (. K2z2  +  F)  (34) 

where  K2  £  R3x3  is  a  positive  constant  diagonal  gain 
matrix,  F(t)  £  R  is  a  neural  network  based  feedforward 
term  which  is  specified  in  a  manner  to  be  bounded1,  and  the 
diagonal  matrix  B  £  R3x3  is  explicitly  given  as 

'  l/6n  0  0 

B  =  0  I/622  0  (35) 

0  0  1/633  . 

where  611,622,633  £  R  are  the  diagonal  elements  of  the 
thrust  configuration  matrix  of  (5).  After  substituting  the 
thrust  input  designed  in  (34)  into  the  open  loop  dynamics 
of  (32),  the  following  closed  loop  error  system  for  z2(t)  is 
developed 

Mz2  =  —K2z2  —  -Mz2  —  zi  +  Ni  +  Nld 

0  612/622  613/633 

+  0  0  623/633  K2z2  (36) 

0  0  0 

where  6y  represents  the  element  from  the  ith  row  and  jth 
column  of  the  matrix  B ,  Nx  =  Ni  —  Nid,  and  IVi(-)  £  R3 
are  defined  as  follows 

N1  =  N-  BBF.  (37) 

1  For  details  of  the  structure  of  the  bounded  neural  network  feedforward 
term,  see  Appendix  II. 


bounded  in  the  following  manner 


In  (36),  the  term  Nid(-)  G  R3  is  defined  as  follows 

Nu  =  (38) 

=  -  BBF.  (39) 

Note  that  from  Property  1,  Nid(t)  can  be  shown  to  be 
bounded  if  u(t),  F(t),  rjd(t),  r}d(t),  rid(t)  G  C (x.  Using  the 
fact  that  the  desired  trajectory  and  its  corresponding  time 
derivatives  are  bounded,  it  can  be  shown  that  N\  can  be 
upper  bounded  in  the  following  manner 

Ni  </5(INI)IN|  (40) 

where  z  =  [zf,zJ]T  G  R6  and  p(-)  is  a  globally  invertible, 
positive,  non-decreasing  function. 

Remark  5:  The  thrust  control  input  defined  in  (34)  re¬ 
quires  knowledge  of  the  diagonal  elements  of  the  thrust 
configuration  matrix.  If  we  expand  (34),  it  can  be  ob¬ 
served  that  the  first  element  of  U(t)  is  defined  as  u\  = 
(1/fen)  (k2iZ2i  +  j\  j ,  where  k’2i  is  the  first  diagonal  ele¬ 
ment  of  A'2,  221(f)  is  the  first  element  of  22(f),  and  /j(f)  is 
the  first  element  of  F(t).  From  the  commutation  strategy  of 
(7a)-(7b),  u\(t)  is  used  solely  to  compute  the  thrust  inputs 
for  the  tugboats  uia(t)  and  u±b(t).  Under  the  assumption 
that  each  tug  knows  its  own  position  along  the  vessel’s 
hull  ( i.e .  bn  will  be  known  by  the  controller  and 

from  the  decomposition  of  the  thrust  configuration  matrix 
presented  in  Section  II,  it  is  easy  to  observe  that  v  \  (f)  can 
be  calculated  without  knowing  any  other  elements  of  the 
thrust  configuration  matrix  and  as  a  result  does  not  require 
information  of  the  other  tugboats  locations. 


B.  Stability  Analysis 

Theorem  2:  For  the  vessel  model  of  (8)-(9),  the  thrust 
control  input  of  (34)  ensures  that  the  position  tracking  error  is 
globally  uniformly  ultimately  bounded  (GUUB)  in  the  sense 
that 

IM*)II>  IMi)||  <  \/A)exp(— A3 t)  +  P1  (41) 

where  A3  G  R+  denote  constants. 

Proof:  In  order  to  illustrate  the  result  of  (41),  the 
following  non-negative  scalar  function  V2  (f)  G  R  is  defined 
as 

V2  =  \#zi  +  i ZzM(v)z2  (42) 

After  taking  the  time  derivative  of  (42),  substituting  from 
(30)  and  (36),  the  following  expression  can  be  obtained  for 
the  time  derivative  of  U2  (t) 

V2  =  -z[ K1z1  -  z% K2z2  +  Z2N1  +  Z2  Nu 

fbi2\  _  fb13\ 

+K22  (  —  I  221^22  +  K23  I  —  I  £21^23 

+k23  Z22Z23  (43) 

where  22,  is  the  ith  element  of  22,  fc2,;  is  the  ith  diagonal 
element  of  A’2  •  The  expression  for  U2  (f)  in  (43)  can  be  upper 
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where  K„  G 


3x3 


is  a  constant  positive  diagonal  gain 


matrix  which  is  selected  such  that  K2  =  5Kn,  kni  is  the  ith 
diagonal  element  of  Kn,  and  Amln  represents  the  minimum 
eigenvalue.  After  completing  the  squares  on  the  bracketed 
terms  of  (44),  the  following  upper  bound  for  the  expression 
in  (44)  can  be  developed 

V2  <  -z[K!Zi  -  knl  |22i|2  -  kn2  I222I2  -  kn3  \z23 12 


^’22  / 

'  bi2 

knl  V 

,&22 

^23  / 

'b23 

kn2  V 

X33 

1 

1 

*22 


*23 


k2 

^23 

knl 


bl3 

b33 


*23 


where  eo  > 

re 

Ui  5^  Amm 


PdMDINr  +  eo  (45) 
||iVld||2  G  M+.  The  expression  of  (45) 


min  {Kn\ 

can  be  further  upper  bounded  in  the  following  manner 


{ATi}||2i||  -  knl  1^21 12  -  kp  |222|2 


—  kq  1 223  |  +  I 


'  +  £0 


(46) 


where  kp ,  kq  are  positive  scalar  constants  which  are  defined 
as  follows 

kp  =  kn 2  —  (k^/knl)  (^12/^22)  , 


kq  =  km  -  {kl3/knl)  (613/^33)  ,  (47) 


km  =  kn3  {^3/ knf)  (^23/^33)  • 

The  control  gains  are  selected  as  follows  to  ensure  that  the 
variables  of  (47)  are  positive 
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(48) 


The  combination  of  the  four  terms  of  (46)  allows  the 
expression  for  V2  (t)  to  be  upper  bounded  as  follows 


1>2  <  —  [Ar  —  p(||z||)]  ||2||2  +  e0 


(49) 


where  z(t)  was  previously  defined  in  (40)  and  Ai  G  M+ 
represents  a  constant  defined  as  follows 


Ai  =  min{ATOj„{Afi},  min{fc„i,  kp,  kq}}. 


(50) 


Provided  that  the  gain  values  are  selected  to  ensure  that 
Ai  >  /^(  ||  2 1| ) ,  then  the  time  derivative  of  V2(i)  can  be  upper 
bounded  as 


V2  <  —  A2H2H2  +  eo 


(51) 


where  A2  =  Ai  —  p(||z||)  £  R+.  From  (42),  (51)  can  be 
upper  bounded  as  follows 

V2  <  — A3V2  +  £0  (52) 


where  A3  £  R+  is  defined  as  follows 

2A2 


A3  — 


(53) 


max{  1 ,  A max{M{v)}}  ' 

The  differential  inequality  of  (52)  can  now  be  utilized  to 
formulate  an  upper  bound  on  V2(t)  as  follows 

v2 (t)  <  V2(0)exp(-X 3t)  +  [1  -  exp(— A3£)]  .(54) 

A3 

Based  on  (42),  (52),  and  the  definition  of  z(t),  (54)  can  be 
rewritten  in  the  following  form 

lk(0H2  <  ftexp(-A3i)  +  (31  (55) 


where  /30,/31  £  R+  are  defined  as 


ft 


max{  1,A  max{M{v)}}  2 

min{l,  \min{M{v)}}  ~ 

_ 2eo _ 

min{  1 ,  \min{M(v)}}\3  ' 


From  (55)  and  the  definition  of  z(t),  the  result  given  in 
(41)  can  be  observed.  Since  z\(t)  £  Coo,  it  is  known  that 
(71  —  V d)  £  £°o’  and  hence,  the  position  tracking  result 
is  achieved.  Given  (41),  (30)  and  (36)  can  be  utilized  to 
prove  that  ii(f),i2(i)  £  Coo ■  After  using  the  fact  that 
z2(t),F(t. )  £  Coo,  the  expression  of  (34)  can  be  employed 
to  prove  that  U(t)  £  Coo ■  ■ 


VI.  Simulation  Results 

In  order  to  verify  operation  and  to  evaluate  performance 
of  the  both  the  adaptive  and  robust  control  strategies,  simula¬ 
tions  were  performed  within  the  MATLAB  SIMULINK  en¬ 
vironment  utilizing  the  model  parameter  values  of  Cybership 
II  [3]  for  the  disabled  vessel  dynamics.  In  order  to  account 
for  the  unidirectional  thrust  input  constraint,  six  tugboats 
were  positioned  at  the  following  locations  with  respect  to 
the  vessel’s  center  of  mass 

Lla  =  (-0.6275, 0),  Llb  =  (0.6275, 0), 

L2a  =  (-0.25,  -0.145),  L2b  =  (0.25, 0.145),  (m) 

L3a  =  (-0.25, 0.145),  L3b  =  (0.25,  -0.145). 

In  addition,  the  incident  angle  of  each  tugboat  with  respect 
to  the  hull  was  selected  according  to  (6)  in  order  to  ensure 
that  the  resulting  configuration  matrix  is  upper  triangular  and 
are  given  as  follows 

Otla  =  0,  Qfih  =  7T,  Oi2a  —  7t/4,  ,  . 

_  ,,  (rad) 

a  26  =  a2a  +  7r,  a3a  =  77T/4,  a3b  =  a3a  —  n. 

The  desired  position  trajectory  for  both  control  strategies 
was  selected  as  follows 

V d  =  [0, 0,  (tt/2)  (1  -  exp(0.1  t))}T 

with  the  initial  position/orientation  of  the  vessel  set  to 
77(0)  =  [1.0, 1.0, 7r/8]T.  The  control  gain  values  were  se¬ 
lected  according  to  the  conditions  set  forth  in  the  stability 


analysis  through  a  trial  and  error  procedure  that  produced  the 
best  overall  position  tracking  performance.  For  the  adaptive 
control  strategy  of  Section  IV,  the  control  gains  were  selected 
in  the  following  manner 

K  =  diag{  20.0, 20.0, 20.0},  a  =  diag{  1.0, 1.0, 1.0}, 

r1  =  20.0I12,  r2  =  2.0-I7 

where  Ii2  £  R12x12  and  I7  £  R7x7  represent  identity 
matrices.  The  initial  estimates  0i(t)  £  S}3,92(t)  £  R7) 
were  initialized  as  follows 

01(0)=#!  ft(0)  =  0.0 

where  0^  was  assumed  to  be  25%  of  the  true  value  6\ .  For  the 
robust  control  strategy,  the  control  gain  values  were  selected 
as  follows 

I<i  =  diag{  1.0, 1.0, 1.0},  K2  =  diag{  500.0, 100.0, 2.0}. 

Figures  2-4  and  Figures  7-9  illustrate  the  position  tracking 
results  obtained  from  adaptive  and  robust  control  strategies, 
respectively.  As  seen  from  the  position  tracking  error  Figures 
2  and  7,  both  controllers  are  able  to  maneuver  the  vessel 
along  the  desired  vessel  tracking.  It  is  also  interesting  to 
observe  the  unidirectional  tugboat  thrust  inputs  of  Figures  4 
and  9.  As  seen  in  these  figures  even  when  the  tracking  error 
signal  is  close  to  zero,  a  small  amount  of  thrust  is  maintained 
to  keep  the  tugboat  in  contact  with  the  vessel  hull  when  not 
in  use  (this  is  the  result  of  the  small  constant  70  present  in  the 
commutation  strategy  of  (7a)  and  (7b).  For  both  simulations 
the  constant  was  selected  as  j0  =  \/5). 

VII.  Conclusion 

In  this  paper,  two  control  strategies  have  been  presented 
for  the  positioning  of  an  unactuated  vessel  by  multiple,  au¬ 
tonomous  tugboats.  Specifically,  an  adaptive  position  control 
strategy  was  developed  to  target  the  problem  where  the  loca¬ 
tion  of  the  tugboats  about  the  vessel  hull  is  not  known.  As  a 
result  of  this  uncertainty,  communication  between  members 
of  the  tugboat  group  is  required  to  adequately  compensate 
for  the  unknown  thrust  matrix.  To  offer  attractiveness  for 
recovery  of  disabled  vessels  within  hostile  environments,  a 
second  robust  control  strategy  is  developed  that  does  not 
require  the  passing  of  data  between  the  tugboats.  To  facil¬ 
itate  this  lack  of  communication/coordination,  the  tugboat 
locations  about  the  vessel  are  selected  a  priori  in  a  manner 
such  that  the  resulting  thrust  configuration  matrix  is  upper 
triangular.  This  structure  allows  for  an  hierarchical  control 
design  approach  where  “upper  level”  tugboats  treat  thrust 
inputs  from  “lower  level”  tugs  as  disturbances.  Both  control 
approaches  have  been  verified  via  simulation.  There  is  much 
work  to  be  considered  for  future  work.  One  of  the  major 
assumptions  of  this  work  is  that  the  tugboat  locations  and 
their  incident  angles  to  the  hull  remain  fixed  throughout  the 
maneuver;  therefore,  subsequent  work  will  include  tugboat 
dynamics  within  the  scope  of  the  control  design.  Also,  the 
control  design  must  also  account  for  the  possibility  of  the 
tugboat  slipping  along  the  hull.  As  a  result,  future  control 
strategies  should  include  the  constraint  that  lateral  forces 


along  the  hull  from  the  tug  remain  below  the  contact  static 
friction  coefficient  [17].  In  addition,  work  is  under  way  for 
the  validation  of  both  control  strategies  via  experiments. 
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Appendix  I 
Proof  of  theorem  I 


defined  as  follows 

V\  =  —rTM*r  +  -61 ri  16i  +  —92  U  1@2-  (56) 

After  taking  the  time  derivative  of  (56),  substituting  from 
(27),  the  following  expression  for  the  time  derivative  of  V\  (t) 
can  be  obtained 

Vi  =  - rTKr  +  |Y1Tr  -  (57) 

+02  [Y?r  -  T2lb2 

After  substituting  the  projection  based  update  laws  of  (21) 
and  (26)  into  (57),  the  dynamics  for  V(t)  can  be  upper 
bounded  in  the  following  manner 

Vi  <  -Xmin{K}\\r\\  (58) 

where  A min{K}  represents  the  minimum  eigenvalue  of  the 
gain  matrix  K. 

From  the  above  expressions  of  (56)  and  (58),  the  signals 
r(f),  9i(t),  92(t)  can  be  observed  to  be  bounded.  Since 
r(t)  G  £oo,  the  definition  of  r(t)  of  (18)  can  then  be  utilized 
to  show  that  position/velocity  tracking  error  signals  are  also 
bounded  (e(t),e(t)  G  Coo).  Since  the  desired  trajectory  was 
selected  such  that  r]d(t),  rjd(t)  G  Coo,  v(t),v(t)  G  Coo- 
Also  since  9i(t),92(t)  G  Coo,  the  parameter  estimates 
9 i{t)  and  92(t)  are  bounded  signals  (note  that  9i(t)  is 
bounded  independent  of  the  stability  analysis  since  the  pro¬ 
jection  based  update  algorithm  of  (23)  was  utilized).  Since 
r(t),e(t)  G  Coo,  Property  1  and  (20)  is  utilized  to  show  that 
Y2(t)  G  Coo',  hence  U(t)  can  be  observed  to  be  a  bounded 
signal  from  (25).  At  this  point,  (21),  (22),  (26),  the  projection 
based  update  algorithm  of  (23),  and  previous  arguments 
can  be  utilized  to  show  that  to  show  that  9i,92  G  Coo- 
Since  e(t),e(t)  G  C2  fl  Coo,  Barbalat’s  Lemma  [18]  can 
be  employed  to  illustrate  the  position  tracking  result  in  the 
sense  that  e(t)  — »  0  as  t  — +  oo.  ■ 

Appendix  II 

Neural  Network  Feedforward 
The  neural  network  feedforward  component  is  computed 
using  a  two  layer  network  using  ‘n’  neurons.  This  feedfor¬ 
ward  term  must  always  be  bounded  (i.e.  /(f)  G  Coo),  to  this 
end  the  neural  network  component  is  defined  as  follows  [14] 

/  =  WTa  (vTx^j  .  (59) 

where  W(t)  G  Rrax3  and  V(t)  G  R13xn  are  es¬ 
timated  weight  matrices,  and  x(t)  G  R13  is  the 
input  vector  to  the  neural  network  which  is  se¬ 
lected  as  i  =  [  1,  r%,  f]d,  rfi,  'rjj  ]  where 

i]d(t),  Vd(t)i  Vdh)  represent  the  desired  trajectory 

and  its  first  three  derivatives.  The  vector  activation 
function  of)  G  R"  i— >  R™  is  defined  as  follows 
a(u>)  =  [a(ixi),  cr(u) 2),  •••  ,  a(ujn))T  G  R"  where  u>  = 
[to i,  uj2,  ■  ■  ■  ,  u}n]T  G  R”  and  o(s)  :  R  i— >  R  is  the  sigmoid 
activation  function  defined  as 


Proof:  In  order  to  illustrate  the  position  tracking  result, 
the  following  non-negative  scalar  function  Vj  (t)  G  R  is 


cr(s) 


1 

1  +  exp(-s) 


(60) 


The  gradient  of  the  vector  activation  function,  denoted  by 
a  (•)  £  R”x"  can  be  expressed  in  closed  form  as  follows, 
[19] 

<7(10)  =  diag{cr{io )}  [I  —  diag{i r(w)}] .  (61) 

To  ensure  that  the  weights  generated  from  the  weight  update 
laws  are  bounded,  the  update  laws  were  defined  as  follows 
[14] 


W  =  —awW  +  7 xa  (vTx^  sat (e2  +  C)T  (62) 


V  =  - avV  +  72x 


f  (vTx^  iTsat  (e2  +  0 


(63) 


where  av,  aw  £  M+  are  small  constants,  7i,  72  G  R+  are 
control  gains  which  effect  the  learning  speed,  the  function 
saf(^)  :  R™  1— >  R™  is  a  saturation  function,  and  the  auxiliary 
signal  £(f)  £  R3  is  a  surrogate  (i.e.  a  dirty  derivative 
operation)  for  the  signal  <22  (t)  which  is  defined  as  follows 


C  =  ~(e2  ~  f])  (64) 

e 

where  e  £  R+  is  a  small  constant,  and  the  signal  77(f)  £  R3 
is  updated  according  to  the  following  expression 

V  =  -(e2-w)-  (65) 

£ 

From  equations  (59)-(65)  and  the  fact  that  the  input  vector 
to  the  neural  network  is  bounded,  it  is  easy  to  show  that  the 
weight  matrices  W (f )  and  V  ( t )  are  bounded,  and  hence,  the 
output  from  the  neural  network,  /(f),  is  bounded. 


Appendix  III 
Simulation  Figures 
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Fig.  3:  Tracking  error  for  the  case  with  full  communication. 
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Fig.  4:  Control  input  for  the 
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Fig.  2:  Actual  (solid  line)  and  desired  (dashed  line)  trajecto¬ 
ries  for  the  case  with  full  communication. 
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Fig.  5:  Dynamic  parameter  estimates  Qi(t),  for  the  case  with 
full  communication. 
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Fig.  8:  Tracking  error  for  the  case  with  no  communication. 
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Fig.  6:  Parameter  estimates  62  (t),  for  the  case  with  full 
communication. 
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Fig.  9:  Control  input  for  the  case  with  no  communication, 
Uia  (solid  line),  Uib  (dashed  line). 


Fig.  7:  Actual  (solid  line)  and  desired  (dashed  line)  trajecto¬ 
ries  for  the  case  with  no  communication. 


